GBSV Mini-Challenge 1 Notebook 2

Day 4

Daten

Für die Analyse benötige ich Daten, welche ein wiederkehrendes Muster aufweisen. In den letzten Tagen habe ich mich mit den Wasserstandsdaten des St.-Lawrence Stroms in Montreal beschäftigt. Diese Daten sind stündlich über mehrere Jahre verfügbar und weisen ein Muster auf, welches Jährlich wiederkehrt.

Alter Hafen von Montreal im Winter

Autokorrelation über alle Jahre

import pandas as pd
import numpy as np
import matplotlib.pyplot as plt
import statsmodels.api as sm
from scipy.signal import correlate, find_peaks
import matplotlib.dates as mdates
import calendar


# Einlesen der CSV-Datei
water_raw = pd.read_csv(
    "data/st_lawrence_water.csv",
    skiprows=8, 
    usecols=[0, 1],                # nur die ersten beiden Spalten laden
    names=["Observed_date", "Sealevel"],    # Spaltennamen setzen
    parse_dates=["Observed_date"],      # Datum parsen
    encoding="latin1"
)

print(water_raw.head())

        Observed_date  Sealevel
0 2015-01-01 00:00:00      0.32
1 2015-01-01 01:00:00      0.32
2 2015-01-01 02:00:00      0.33
3 2015-01-01 03:00:00      0.33
4 2015-01-01 04:00:00      0.33

# Plotten der Daten einmal über die ganzen Jahre und für ein Jahr
plt.figure(figsize=(16, 10))
plt.subplot(2, 1, 1)
plt.plot(water_raw["Observed_date"], water_raw["Sealevel"])
plt.title("Wasserstand über die Jahre")
plt.xlabel("Datum")
plt.ylabel("Wasserstand (m)")
plt.grid()
plt.xlim([water_raw["Observed_date"].min(), water_raw["Observed_date"].max()])
plt.ylim([water_raw["Sealevel"].min(), water_raw["Sealevel"].max()])
plt.subplot(2, 1, 2)
plt.plot(water_raw["Observed_date"], water_raw["Sealevel"])
plt.title("Wasserstand im Jahr 2016")
plt.xlabel("Datum")
plt.ylabel("Wasserstand (m)")
plt.grid()
plt.xlim([pd.Timestamp("2016-01-01"), pd.Timestamp("2017-01-01")])
plt.ylim([water_raw["Sealevel"].min(), water_raw["Sealevel"].max()])

(-0.29, 3.29)

Plot Optimiert

# Anzahl der Lags
lags = len(water_raw["Sealevel"])-1

# Berechnung und Plot der Autokorrelationsfunktion
sm.graphics.tsa.plot_acf(water_raw["Sealevel"], lags=lags, title="Correlogram - Stündliche Daten")
plt.xlabel("Lags in Stunden")
plt.ylabel("Autokorrelationskoeffizient")
plt.show()

Um das ganze etwas klarer zu machen und auch mehr die saisonalität zu sehen, habe ich die Daten auf monatliche Intervalle gemittelt.

water_monthly = water_raw.resample('ME', on='Observed_date').mean()

Plot optimiert

# Autokorrelationsfunktion für monatliche Daten
lags = len(water_monthly)-1
sm.graphics.tsa.plot_acf(water_monthly["Sealevel"], lags=lags, title="Correlogram - Monatlich Gemittelt")
plt.xlabel("Lags in Monaten")
plt.ylabel("Autokorrelationskoeffizient")
plt.show()

# Periodizität der Daten aufzeigen

# Autokorrelation für monatliche Daten
acf_values = sm.tsa.stattools.acf(water_monthly["Sealevel"], nlags=lags)
lags_array = np.arange(len(acf_values))

# Peaks finden (ausser bei Lag = 0)
peaks, _ = find_peaks(acf_values, height=0.2)

plt.figure(figsize=(10,4))
plt.plot(lags_array, acf_values, label='Autokorrelation')
plt.scatter(lags_array[peaks], acf_values[peaks], color='red', label='Entdeckte Peaks')
plt.title("Autokorrelation Peaks - Monatliche Daten")
plt.xlabel("Lag (Monate)")
plt.ylabel("Autokorrelation Koefficient")
plt.grid(True)
plt.legend()
plt.show()

# Periodendauer in Monaten
if len(peaks) > 1:
    period_months = np.mean(np.diff(peaks))
    print(f"Geschätzte Periodendauer: {period_months:.1f} Monate (~{period_months/12:.2f} Jahre)")

Geschätzte Periodendauer: 20.0 Monate (~1.67 Jahre)

Im Auto-Korrelationsdiagramm sieht man nun, dass nach ca. 12 Monaten die Korrelation wieder sehr hoch ist, was auf eine jährliche Saisonalität hinweist. Dies wird auch erkenntlich, wenn man nach den höchsten Punkten (Peaks) im Diagramm sucht. Die Periodendauer wird dabei auf 12 errechnet, was durchaus Sinn ergibt. Um das Ganze noch etwas einfacher darzustellen, kann mach auch die Daten in 4 Quartale pro Jahr aufteilen. Dies schwächt natürlich die Genauigkeit etwas ab, macht das Muster aber noch klarer ersichtlich.

Plot optimiert

water_quarterly = water_raw.resample('QE', on='Observed_date').mean()

# Autokorrelationsfunktion für monatliche Daten
lags = len(water_quarterly)-1
sm.graphics.tsa.plot_acf(water_quarterly["Sealevel"], lags=lags, title="Correlogram - Quartalsmässig Gemittelt")
plt.xlabel("Lags in Quartalen")
plt.show()

Hier sieht man noch genauer, jedes Vierte Quartal weisst wieder eine höhere Korrelation auf. Nach 4 Jahren ist die Korrelation in allen Plots nahe bei 0. Wenn die jeweils 4 Aufeinanderfolgende Jahre betrachtet werden, sieht man, dass diese oft unterschiedlich zu anderen 4 Jahreszeiträumen sind. Wahrscheinlich hebt sich dieser Effekt im Anschluss bei 5 und 6 Jährigen Zeiträumen wieder etwas auf.

Day 5

Problemstellung

Die wiederholenden Muster in den Daten sind wichtig zu verstehen. Es ist enorm wichtig, diese Muster zu erkennen, um Vorhersagen über zukünftige Ereignisse treffen zu können. In diesem Fall muss bei der Planung von Bauprojekten wie Brücken oder Hafenanlagen der Wasserstand berücksichtigt werden. Und dies sollte für alle Jahreszeiten und nicht nur für den Durchschnittswasserstand über das Jahr hinweg geschehen. Ansonsten könnten die Brücken im Frühling/Sommer zu niedrig sein oder die Hafenanlagen im Herbst/Winter zu hoch. Zudem ist wichtig zu erkennen, dass die Wasserstände von Jahr zu Jahr über mehrere Jahre hinweg zwar ähnlich sind, aber nicht Identisch. Beispielsweise ist bei Vier Jahren in folge fast keine Autokorrelation mehr vorhanden. Es muss also beim Bau auch darauf geachtet werden, dass die Wasserstände von Jahr zu Jahr etwas variieren können.

Day 6

Observation

Die Autokorrelationsanalyse zeigt deutlich, dass nach etwa 12 Monaten eine starke Korrelation besteht. Das bestätigt den saisonalen Jahreszyklus im Wasserstand, der für Hafen- und Brückenplanung besonders relevant ist. Gleichzeitig fällt auf, dass sich die Autokorrelation nach 4 Jahren fast auf null reduziert, was auf natürliche Variabilität zwischen den Jahren hinweist.

Interpretation

Für den Use Case bedeutet dies, dass Planungen nicht nur auf einem Jahresmittel basieren dürfen. Konstruktionen müssen Schwankungen zwischen verschiedenen Jahren berücksichtigen, da die Wiederholbarkeit über längere Zeiträume abnimmt. Eine reine Orientierung am Durchschnittswert könnte zu Fehlplanungen führen, z. B. zu niedrig bemessenen Brücken im Frühjahr.

Diskussion

Die Ergebnisse sind im Kontext des Anwendungsfalls sehr wertvoll: Einerseits kann die jährliche Saisonalität klar für Prognosen genutzt werden, andererseits zeigen die niedrigen Korrelationen nach mehreren Jahren, dass die Variabilität zusätzliche Sicherheitsreserven in der Bauplanung notwendig macht. Kritisch zu hinterfragen ist, ob externe Einflüsse (z. B. Klimawandel, Niederschläge) die Korrelation langfristig weiter schwächen dies ist ein Aspekt, der in zukünftigen Analysen berücksichtigt werden muss.

Methodenwahl

Ich habe die Autokorrelationsfunktion (ACF) gewählt, da sie eine etablierte Methode ist, um periodische Muster in Zeitreihen zu erkennen. Sie ist besonders geeignet für saisonale Daten, wie Wasserstände, da sie Abhängigkeiten über Monate und Jahre hinweg direkt sichtbar macht. Alternativen wie Fourier-Analysen wären ebenfalls möglich gewesen, fokussieren aber stärker auf Frequenzen, während ACF die Zeitabhängigkeit klarer darstellt.

Parameterwahl

Ich habe die Daten auf monatliche und quartalsweise Intervalle gemittelt. Dadurch werden kurzfristige Schwankungen (z. B. Wetter oder Gezeiten) geglättet und der saisonale Zyklus tritt klarer hervor. Lags bis zur Länge der Zeitreihe wurden berücksichtigt, um auch mehrjährige Muster sichtbar zu machen. Die Quartalsbetrachtung wurde gewählt, weil sie einen guten Kompromiss zwischen Glättung und Erhalt saisonaler Details darstellt. Diese Quartalsbetrachtung ist im Rückblick auf das vorhergehende Notebook und die dortige Erkenntniss, dass eine Abtastrate von 3 Monaten Sinnvoll ist, entstanden.

Min und Max Berechnen

min_level = water_raw["Sealevel"].min()
max_level = water_raw["Sealevel"].max()
amplitude = max_level - min_level

print(f"Minimaler Wasserstand: {min_level:.2f} m")
print(f"Maximaler Wasserstand: {max_level:.2f} m")
print(f"Jahres-Amplitude: {amplitude:.2f} m")

water_raw["Year"] = water_raw["Observed_date"].dt.year
yearly_range = water_raw.groupby("Year")["Sealevel"].agg(lambda x: x.max() - x.min())
print(f"Mittlere jährliche Amplitude: {yearly_range.mean():.2f} m (±{yearly_range.std():.2f})")

Minimaler Wasserstand: -0.29 m
Maximaler Wasserstand: 3.29 m
Jahres-Amplitude: 3.58 m
Mittlere jährliche Amplitude: 2.54 m (±0.43)

Min und Max

Die Minima und die Maxima zeigen klar, das sich der Wasserstand veränder und zwar enorm stark. Es ist ein Unterschied von 3.58 m zwischen dem höchsten und dem niedrigsten Wasserstand. Dies ist enorm viel und muss bei der Planung von Brücken und Hafenanlagen unbedingt berücksichtigt werden. Mithilfe der mittleren Jahresamplitude von 2.54 m wird zudem klar, dass die Wasserstände im Jahresverlauf stark schwanken. Gerade für die etwas niedrigeren Brücken wie beispielsweise die Pont de la Concorde ist dies für die Stützten wichtig zu wissen, dass der höchste Wasserstand fast 4 Meter höher ist als der niedrigste. Für die ganz grossen Brücken wie die Pont Jacques-Cartier ist dies zwar auch wichtig, aber da die Brücke sehr hoch ist, ist der Unterschied von 3.58 m im Verhältnis zur Gesamthöhe der Brücke nicht mehr so relevant.

Day 7

Use Case

Im Frühjahr steigt der Wasserstand im St.-Lawrence-Strom stark an, weil Schneeschmelze und erhöhte Zuflüsse einsetzen. Dieser Anstieg kann zu Einschränkungen der Schifffahrt, Überflutungen in Uferzonen und einer Belastung der Hafeninfrastruktur führen. Für die Hafenplanung und den Hochwasserschutz ist es daher entscheidend, den Zeitpunkt des beginnenden Wasseranstiegs frühzeitig zu erkennen und mit früheren Jahren zu vergleichen.

Auswahl des Snippets

Als Referenz habe ich den Zeitraum vom 15. Februar bis 15. Mai 2017 gewählt. Dieser Zeitraum umfasst den typischen Anstieg des Wasserstands im Frühling und den Beginn des Sommers. Im Jahr 2017 war der Anstieg besonders ausgeprägt, dies ist gut, da wir vor allem sehr hohe Wasserstände erkennen wollen. Bei den anderen Jahren wird dann jeweils der Zeitraum zwischen Januar und Juni untersucht, um ähnliche Muster zu finden. Über das gesamte Signal wird eine gleitende Pearson-Korrelation mit dem Snippet berechnet. So können ähnliche Wasserstandsmuster in anderen Jahren erkannt werden. Für jedes Jahr wird der Zeitpunkt mit der höchsten Korrelation ausgewählt, um die besten Übereinstimmungen zu identifizieren.

Monate beim Plot optimiert und neue Kennzahlen für Monate

# prepare daily signal
sig = water_raw.set_index('Observed_date')['Sealevel'].asfreq('D').interpolate(limit=16)

# 1) choose snippet (Referenzjahr)
snippet = sig.loc['2017-02-15':'2017-05-15'].copy()

# 2) prepare data
x = sig.values.astype(float)
t = snippet.values.astype(float)
L = len(t)

# 3) compute sliding Pearson
dot_valid = correlate(x, t, mode='valid')  # length N-L+1
s = pd.Series(x, index=sig.index)
mean_x = s.rolling(window=L).mean().to_numpy()[L-1:]
std_x  = s.rolling(window=L).std(ddof=0).to_numpy()[L-1:]
mean_t = np.mean(t)
std_t = np.std(t)
eps = 1e-8
pearson = (dot_valid - L * mean_x * mean_t) / (L * (std_x * std_t + eps))

timestamps_valid = sig.index[:len(pearson)]

# 4) pick top match per year in Jan-Jun window
matches = []
for y in range(sig.index.year.min(), sig.index.year.max()+1):
    mask_year = (timestamps_valid.year == y) & (timestamps_valid.month >= 1) & (timestamps_valid.month <= 6)
    if mask_year.any():
        idxs = np.where(mask_year)[0]
        best_idx = idxs[np.argmax(pearson[idxs])]
        matches.append((timestamps_valid[best_idx], pearson[best_idx]))

matches = [m for m in matches if m[1] > 0.2]  # optional absolute filter

print("Yearly best matches (Jan-Jun):")
for t, c in matches:
    print(t.date(), round(c,3))

# 5) plot pearson with yearly markers (alle 2 Monate + 90° Rotation)
plt.figure(figsize=(12,4))
plt.plot(timestamps_valid, pearson, label='Pearson-Korrelation', color='tab:blue')
for t, c in matches:
    plt.scatter(t, c, color='red', zorder=5)

plt.title('Sliding Pearson – Übereinstimmungen in Jan–Jun')
plt.xlabel('Monat')
plt.ylabel('Pearson-Korrelation')
plt.grid(True)

# Monatslabels alle 2 Monate
plt.gca().xaxis.set_major_locator(mdates.MonthLocator(bymonth=[1,3,5,7,9,11]))
plt.gca().xaxis.set_major_formatter(mdates.DateFormatter('%b %Y'))
plt.xticks(rotation=90)
plt.tight_layout()
plt.show()

# 6) KPI: Häufigkeit pro Monat der besten Übereinstimmung
months = pd.Series([m[0].month for m in matches])
month_counts = months.value_counts().sort_index()

print("\nHäufigkeit des Startmonats der besten Korrelation pro Jahr:")
for month_num, count in month_counts.items():
    month_name = calendar.month_name[month_num]
    print(f"  {month_name:<10}: {count} Jahr(e)")

# 7) overlay example: show match overlays for a few years (alle 2 Monate)
plt.figure(figsize=(12,6))
plt.plot(sig.index, sig.values, alpha=0.4, label='Gesamtsignal', color='gray')
for t, c in matches:
    start = sig.index.get_loc(t)
    end = start + L
    plt.plot(sig.index[start:end], sig.values[start:end], linewidth=2, label=f'{t.year} (corr={c:.2f})')

plt.legend(loc='center left', bbox_to_anchor=(1, 0.5))
plt.title('Signal mit überlagerten besten Jahres-Übereinstimmungen')
plt.xlabel('Monat')
plt.ylabel('Wasserstand (m)')
plt.grid(True)

# Monatslabels alle 2 Monate
plt.gca().xaxis.set_major_locator(mdates.MonthLocator(bymonth=[1,3,5,7,9,11]))
plt.gca().xaxis.set_major_formatter(mdates.DateFormatter('%b %Y'))
plt.xticks(rotation=90)
plt.tight_layout()
plt.show()

Yearly best matches (Jan-Jun):
2015-05-24 0.662
2016-01-01 0.915
2017-02-15 1.0
2018-03-06 0.585
2019-02-24 0.899
2020-01-21 0.879
2021-05-18 0.622
2022-01-29 0.872
2023-02-15 0.878
2024-02-17 0.676

Häufigkeit des Startmonats der besten Korrelation pro Jahr:
  January   : 3 Jahr(e)
  February  : 4 Jahr(e)
  March     : 1 Jahr(e)
  May       : 2 Jahr(e)

Auf dem Plot sieht man gut, dass die Wasseranstiege im Frühling in den verschiedenen Jahren gut erkannt werden. Gerade in den Jahren, in dennen ein starker Anstieg stattgefunden hat, ist die Korrelation sowie die Übereinstimmung mit dem Snippet besonders hoch. Dies macht auch sinn, da dies im Jahr 2017 auch der Fall war. In den Jahren, in denen der Anstieg weniger stark war, ist die Korrelation entsprechend geringer und passt Teilweise schlechter zum Snippet. Da aber gerade hohe Wasserstände relevant sind, ist dies in Ordnung.

Die Schneeschmelze hängt stark von den Lufttemperaturen ab. In Jahren mit warmem Frühling verschiebt sich der Wasseranstieg entsprechend nach vorne. Eine Erweiterung dieser Analyse könnte sein, Temperaturdaten aus Montreal oder dem Oberlauf des St.-Lawrence-Stroms zu berücksichtigen, um den zeitlichen Zusammenhang zwischen Temperaturanstieg und Wasserstand zu quantifizieren. Gerade die Temperaturen bei den Grossen Seen wäre hier interessant, da diese den Zufluss in den St.-Lawrence-Strom stark beeinflussen.

Day 8

Änderung der Amplitude

Die Schneeschmelze und der damit verbundene Wasserstand im St.-Lawrence-Strom können von Jahr zu Jahr variieren. Faktoren wie die Menge des gefallenen Schnees im Winter, die Temperaturentwicklung im Frühling und Niederschlagsmuster beeinflussen die Höhe des Wasserstands. Das Muster bleibt jedoch ähnlich, nur die höhe der Ausschläge ändert sich.

Leichte Zeitliche Verschiebung (Phase shift)

Die Schmelze kann früher oder später einsetzen, abhängig von den Wetterbedingungen. Ein wärmerer Frühling kann zu einem früheren Anstieg des Wasserstands führen, während ein kälterer Frühling den Anstieg verzögern kann. Dies führt zu einer zeitlichen Verschiebung des Musters ohne das grundsätzliche Muster zu verändern.

Plots geupdatet und KPI hinzugefügt

# Original snippet (Referenz)
snippet = sig.loc['2017-02-15':'2017-05-15'].copy()
t = snippet.values.astype(float)
L = len(t)

# 1) Transformationen: amplitude scaling + phase shift
amplitude_factor = 1.3      # 30% höhere Schneeschmelze
shift_days = 14              # 2 Woche später

sig_modified = sig.copy()
sig_modified.iloc[shift_days:] = sig.iloc[:-shift_days].values  # simple shift
sig_modified = sig_modified * amplitude_factor  # scale amplitude

# 2) Compute sliding Pearson cross-correlation
x = sig_modified.values.astype(float)
dot_valid = correlate(x, t, mode='valid')
s = pd.Series(x, index=sig_modified.index)
mean_x = s.rolling(window=L).mean().to_numpy()[L-1:]
std_x  = s.rolling(window=L).std(ddof=0).to_numpy()[L-1:]
mean_t = np.mean(t)
std_t = np.std(t)
eps = 1e-8
pearson = (dot_valid - L * mean_x * mean_t) / (L * (std_x * std_t + eps))
timestamps_valid = sig_modified.index[:len(pearson)]

# 3) Top match per year in Feb-Jun window
matches_modified = []
for y in range(sig_modified.index.year.min(), sig_modified.index.year.max()+1):
    mask_year = (timestamps_valid.year == y) & (timestamps_valid.month >= 2) & (timestamps_valid.month <= 6)
    if mask_year.any():
        idxs = np.where(mask_year)[0]
        best_idx = idxs[np.argmax(pearson[idxs])]
        matches_modified.append((timestamps_valid[best_idx], pearson[best_idx]))

matches_modified = [m for m in matches_modified if m[1] > 0.2]  # optional absolute threshold

# Schönere Print-Ausgabe
print("Yearly best matches (after transformation):")
for t, c in matches_modified:
    print(f"{t.year}-{t.month:02d}: Pearson = {c:.3f}")

# 4) Visualize Pearson with yearly markers (alle 2 Monate, 90°)
plt.figure(figsize=(12,4))
plt.plot(timestamps_valid, pearson, label='Pearson-Korrelation', color='tab:blue')
for t, c in matches_modified:
    plt.scatter(t, c, color='red', zorder=5)

plt.title('Sliding Pearson per year after amplitude & shift')
plt.xlabel('Monat')
plt.ylabel('Pearson-Korrelation')
plt.grid(True)

# Monatslabels alle 2 Monate
plt.gca().xaxis.set_major_locator(mdates.MonthLocator(bymonth=[2,4,6,8,10,12]))
plt.gca().xaxis.set_major_formatter(mdates.DateFormatter('%b %Y'))
plt.xticks(rotation=90)
plt.tight_layout()
plt.show()

# 5) Overlay example: show match overlays for all years
plt.figure(figsize=(12,6))
plt.plot(sig_modified.index, sig_modified.values, alpha=0.4, label='modified signal', color='gray')
for t, c in matches_modified:
    start = sig_modified.index.get_loc(t)
    end = start + L
    plt.plot(sig_modified.index[start:end], sig_modified.values[start:end], linewidth=2, label=f'{t.year}-{t.month:02d} (corr={c:.2f})')

plt.legend(loc='center left', bbox_to_anchor=(1, 0.5))
plt.title('Signal with yearly matches overlaid (after transformation)')
plt.xlabel('Monat')
plt.ylabel('Wasserstand (m)')
plt.grid(True)

# Monatslabels alle 2 Monate
plt.gca().xaxis.set_major_locator(mdates.MonthLocator(bymonth=[2,4,6,8,10,12]))
plt.gca().xaxis.set_major_formatter(mdates.DateFormatter('%b %Y'))
plt.xticks(rotation=90)
plt.tight_layout()
plt.show()

# 6) KPI: Häufigkeit pro Monat der besten Übereinstimmung
months = pd.Series([m[0].month for m in matches_modified])
month_counts = months.value_counts().sort_index()

print("\nHäufigkeit des Startmonats der besten Korrelation pro Jahr:")
for month_num, count in month_counts.items():
    month_name = calendar.month_name[month_num]
    print(f"  {month_name:<10}: {count} Jahr(e)")

Yearly best matches (after transformation):
2015-06: Pearson = 0.662
2016-02: Pearson = 0.752
2017-03: Pearson = 1.000
2018-03: Pearson = 0.585
2019-03: Pearson = 0.899
2020-02: Pearson = 0.879
2021-06: Pearson = 0.622
2022-02: Pearson = 0.872
2023-03: Pearson = 0.878
2024-03: Pearson = 0.676

Häufigkeit des Startmonats der besten Korrelation pro Jahr:
  February  : 3 Jahr(e)
  March     : 5 Jahr(e)
  June      : 2 Jahr(e)

Auch nach den beiden Transformationen (Amplitude und Phase) ist die Korrelation immer noch relativ hoch, was zeigt, dass die grundlegenden Muster im Wasserstand erhalten bleiben. Dies ist wichtig für die Planung, da es bedeutet, dass trotz natürlicher Variabilität ähnliche Muster erkannt und genutzt werden können. Die Korrelationen sind teilweise gleich, teilweise etwas geringer oder höher. Man sieht aber, dass die Muster selbst nicht verändert hat. Die Zeitliche Verschiebung hat bis auf die von mir gesetzten Grenzen keinen einfluss auf die Korrelation und die Amplitudenänderung hat nur einen geringen Einfluss da in der Korrelation die Amplitude ja sowieso normalisiert wird.

Was man aber sieht, ist dass es keine Schmelzen gibt welche im Januar beginnen. Dies wird mit dem Phase shift nach hinten zu tun haben.

Neue Plots zum vergleich zu Original hinzugefügt

# --- NEUER PLOT: Original vs. Modifiziert ---
plt.figure(figsize=(12,5))
plt.plot(sig.index, sig.values, label='Original', alpha=0.6)
plt.plot(sig_modified.index, sig_modified.values, label=f'Modifiziert (×{amplitude_factor}, +{shift_days} Tage)', alpha=0.8)
plt.title('Vergleich: Original vs. Modifiziertes Signal')
plt.xlabel('Monat')
plt.ylabel('Wasserstand (m)')
plt.gca().xaxis.set_major_locator(mdates.MonthLocator(bymonth=[2,4,6,8,10,12]))
plt.gca().xaxis.set_major_formatter(mdates.DateFormatter('%b %Y'))
plt.xticks(rotation=90)
plt.legend()
plt.grid()
plt.tight_layout()
plt.show()

# --- Zusatzplot (Zoom auf Zeitraum des Snippets) ---
plt.figure(figsize=(12,4))
plt.plot(sig.loc['2017-02-01':'2017-06-01'], label='Original', alpha=0.6)
plt.plot(sig_modified.loc['2017-02-01':'2017-06-01'], label='Modifiziert', alpha=0.8)
plt.title('Zoom: Original vs. Modifiziert (Frühjahrsperiode)')
plt.xlabel('Monat')
plt.ylabel('Wasserstand (m)')
plt.gca().xaxis.set_major_locator(mdates.MonthLocator())
plt.gca().xaxis.set_major_formatter(mdates.DateFormatter('%b %Y'))
plt.xticks(rotation=90)
plt.legend()
plt.grid()
plt.tight_layout()
plt.show()

Bei den Vergleichsplots zum Original sieht man gut, dass die Muster trotz der Transformationen sehr ähnlich bleiben. Sie verschieben sich lediglich etwas nach rechts und sind in der Amplitude skaliert. Dies bestätigt, dass die grundlegenden saisonalen Muster im Wasserstand erhalten bleiben, auch wenn sie in ihrer Höhe und ihrem zeitlichen Auftreten variieren. Für die Planung ist dies ein gutes Zeichen, da dies die am natürlichst eintretenden Ereigenisse sein könnten und diese erkannt werden. Anders würde es aussehen, wenn das Muster komplett verändert würde, dafür müsste sich aber das komplette Klima ändern. Dies ist zwar durch den Klimawandel möglich, aber in den nächsten Jahren eher unwahrscheinlich. Um dies noch etwas weiter zu Explorieren werden die Schritte nachfolgend noch einmal unabhängig von einander aber mit verschiedenen Parametern durchgeführt.

Neue Plots mit Parameter Variationen hinzugefügt

# --- Parameterreihe für Amplitudenänderung ---
amplitude_factors = [0.6, 1.0, 1.8, 2.5]
pearson_results_amp = {}

t = snippet.values.astype(float).ravel()
mean_t, std_t = np.mean(t), np.std(t)
eps = 1e-8

plt.figure(figsize=(12,5))

# verschiedene Linienstile zur besseren Unterscheidung
styles = ['-', '--', '-.', ':']

for i, factor in enumerate(amplitude_factors):
    # 1) Skalieren
    sig_amp = sig * factor

    # 2) Sliding Pearson berechnen
    x = sig_amp.values.astype(float).ravel()
    dot_valid = correlate(x, t, mode='valid')
    s = pd.Series(x, index=sig_amp.index)
    mean_x = s.rolling(window=L).mean().to_numpy()[L-1:]
    std_x  = s.rolling(window=L).std(ddof=0).to_numpy()[L-1:]

    pearson = (dot_valid - L * mean_x * mean_t) / (L * (std_x * std_t + eps))
    timestamps_valid = sig_amp.index[:len(pearson)]

    # 3) Speichern für spätere Auswertung
    pearson_results_amp[factor] = {
        "pearson": pearson,
        "timestamps": timestamps_valid
    }

    # Plotten mit unterschiedlichen Linienstilen
    plt.plot(
        timestamps_valid,
        pearson,
        styles[i % len(styles)],
        label=f'Amplitude ×{factor}',
        linewidth=2
    )

plt.title('Sliding Pearson – Effekt verschiedener Amplitudenänderungen')
plt.xlabel('Monat')
plt.ylabel('Pearson-Korrelation')
plt.gca().xaxis.set_major_locator(mdates.MonthLocator(bymonth=[2,4,6,8,10,12]))
plt.gca().xaxis.set_major_formatter(mdates.DateFormatter('%b %Y'))
plt.xticks(rotation=90)
plt.legend(title='Amplitude-Faktor')
plt.grid(True)
plt.tight_layout()
plt.show()

# --- 4) Beispielplot: Original vs. verstärktes Signal ---
plt.figure(figsize=(12,5))
for factor in amplitude_factors:
    plt.plot(sig.index, (sig * factor).values, label=f'×{factor}', alpha=0.7)
plt.title('Originalsignal bei verschiedenen Amplitudenfaktoren')
plt.xlabel('Monat')
plt.ylabel('Wasserstand (m)')
plt.gca().xaxis.set_major_locator(mdates.MonthLocator(bymonth=[2,4,6,8,10,12]))
plt.gca().xaxis.set_major_formatter(mdates.DateFormatter('%b %Y'))
plt.xticks(rotation=90)
plt.legend()
plt.grid(True)
plt.tight_layout()
plt.show()

# --- 5) Quantitative Analyse der Korrelation (Peak-Werte) ---
print("=== Peak-Korrelationen je Amplitudenfaktor ===")
peak_stats = []
for factor, data in pearson_results_amp.items():
    p = data["pearson"]
    max_corr = np.nanmax(p)
    mean_corr = np.nanmean(p)
    median_corr = np.nanmedian(p)
    peak_stats.append((factor, max_corr, mean_corr, median_corr))
    print(f"Faktor ×{factor:>3}:  Max={max_corr:.3f},  Mittel={mean_corr:.3f},  Median={median_corr:.3f}")

# --- 6) Visualisierung der Peak-Korrelationen ---
factors, max_vals, mean_vals, med_vals = zip(*peak_stats)

plt.figure(figsize=(8,5))
plt.plot(factors, max_vals, 'o-', label='Max. Pearson')
plt.plot(factors, mean_vals, 's--', label='Mittelwert Pearson')
plt.plot(factors, med_vals, 'd-.', label='Median Pearson')
plt.title('Änderung der Pearson-Korrelation mit Amplitudenfaktor')
plt.xlabel('Amplitude-Faktor')
plt.ylabel('Korrelationswert')
plt.grid(True)
plt.legend()
plt.tight_layout()
plt.show()

=== Peak-Korrelationen je Amplitudenfaktor ===
Faktor ×0.6:  Max=1.000,  Mittel=0.003,  Median=0.096
Faktor ×1.0:  Max=1.000,  Mittel=0.003,  Median=0.096
Faktor ×1.8:  Max=1.000,  Mittel=0.003,  Median=0.096
Faktor ×2.5:  Max=1.000,  Mittel=0.003,  Median=0.096

Bei der Analyse über die änderung der Amplitude wird klar, dass sich mit der Pearson-Korrelation nichts ändert. Dies ist logisch, da die Pearson-Korrelation die Amplitude ja sowieso normalisiert. Es ist also egal, ob die Amplitude grösser oder kleiner wird, die Korrelation bleibt gleich. Dies zeigt, dass die Mustererkennung robust gegenüber Änderungen in der Amplitude ist. Für die Planung bedeutet dies, dass auch wenn der Wasserstand von Jahr zu Jahr variiert, ähnliche Muster erkannt werden können.

# --- Parameterreihe für Zeitverschiebung ---
shift_days_list = [0, 14, 60, 180]
pearson_results_shift = {}

plt.figure(figsize=(12,5))

for shift_days in shift_days_list:
    # 1) Zeitlich verschieben
    sig_shift = sig.copy()
    if shift_days > 0:
        sig_shift.iloc[shift_days:] = sig.iloc[:-shift_days].values
    elif shift_days < 0:
        sig_shift.iloc[:shift_days] = sig.iloc[-shift_days:].values  # falls rückwärts
    # sonst: shift_days == 0 → Original

    # 2) Sliding Pearson berechnen
    x = sig_shift.values.astype(float).ravel()
    dot_valid = correlate(x, t, mode='valid')
    s = pd.Series(x, index=sig_shift.index)
    mean_x = s.rolling(window=L).mean().to_numpy()[L-1:]
    std_x  = s.rolling(window=L).std(ddof=0).to_numpy()[L-1:]
    mean_t, std_t = np.mean(t), np.std(t)
    eps = 1e-8
    pearson = (dot_valid - L * mean_x * mean_t) / (L * (std_x * std_t + eps))
    timestamps_valid = sig_shift.index[:len(pearson)]

    # 3) Speichern & plotten
    pearson_results_shift[shift_days] = {
        "pearson": pearson,
        "timestamps": timestamps_valid
    }

    plt.plot(timestamps_valid, pearson, label=f'Shift +{shift_days} Tage')

plt.title('Sliding Pearson – Effekt verschiedener Zeitverschiebungen')
plt.xlabel('Monat'); plt.ylabel('Pearson-Korrelation')
plt.gca().xaxis.set_major_locator(mdates.MonthLocator(bymonth=[2,4,6,8,10,12]))
plt.gca().xaxis.set_major_formatter(mdates.DateFormatter('%b %Y'))
plt.xticks(rotation=90)
plt.legend(title='Zeitverschiebung'); plt.grid(True)
plt.tight_layout()
plt.show()

# --- 4) Beispielplot: Original vs. verschobene Signale ---
plt.figure(figsize=(12,5))
for shift_days in shift_days_list:
    sig_shift = sig.copy()
    if shift_days > 0:
        sig_shift.iloc[shift_days:] = sig.iloc[:-shift_days].values
    plt.plot(sig.index, sig_shift.values, label=f'+{shift_days} Tage', alpha=0.7)
plt.title('Originalsignal bei verschiedenen Zeitverschiebungen')
plt.xlabel('Monat'); plt.ylabel('Wasserstand (m)')
plt.gca().xaxis.set_major_locator(mdates.MonthLocator(bymonth=[2,4,6,8,10,12]))
plt.gca().xaxis.set_major_formatter(mdates.DateFormatter('%b %Y'))
plt.xticks(rotation=90)
plt.legend(); plt.grid(True)
plt.tight_layout()
plt.show()

# --- 5) Quantitative Analyse der Korrelation (Peak-Werte) ---
print("=== Peak-Korrelationen je Zeitverschiebung ===")
shift_stats = []
for shift_days, data in pearson_results_shift.items():
    p = data["pearson"]
    max_corr = np.nanmax(p)
    mean_corr = np.nanmean(p)
    median_corr = np.nanmedian(p)
    shift_stats.append((shift_days, max_corr, mean_corr, median_corr))
    print(f"Shift +{shift_days:>3} Tage:  Max={max_corr:.3f},  Mittel={mean_corr:.3f},  Median={median_corr:.3f}")

# --- 6) Visualisierung der Korrelationen vs. Zeitverschiebung ---
shifts, max_vals, mean_vals, med_vals = zip(*shift_stats)

plt.figure(figsize=(8,5))
plt.plot(shifts, max_vals, 'o-', label='Max. Pearson')
plt.plot(shifts, mean_vals, 's--', label='Mittelwert Pearson')
plt.plot(shifts, med_vals, 'd-.', label='Median Pearson')
plt.title('Änderung der Pearson-Korrelation mit Zeitverschiebung')
plt.xlabel('Zeitverschiebung (Tage)')
plt.ylabel('Korrelationswert')
plt.grid(True)
plt.legend()
plt.tight_layout()
plt.show()

=== Peak-Korrelationen je Zeitverschiebung ===
Shift +  0 Tage:  Max=1.000,  Mittel=0.003,  Median=0.096
Shift + 14 Tage:  Max=1.000,  Mittel=0.000,  Median=0.090
Shift + 60 Tage:  Max=1.000,  Mittel=0.007,  Median=0.103
Shift +180 Tage:  Max=1.000,  Mittel=0.018,  Median=0.126

Bei einer sehr starken Änderung der Amplitude (etwa um zwei Monate), sieht man dass sich auch die Korrelation etwas ändert. Dies wird damit zu tun haben, dass wir das Muster gezielt im Frühlich gesucht haben und sich bei einer sehr starken Änderung der Zeitverschiebung das Muster nicht mehr so gut mit dem Snippet deckt. Jedoch wird auch hier klar, dass die Mustererkennung relativ robust gegenüber Zeitlichen Verschiebungen ist. Für die Planung bedeutet dies, dass auch wenn der Wasserstand von Jahr zu Jahr etwas früher oder später ansteigt, ähnliche Muster erkannt werden können.

Day 9

Die Cross-Korrelation hat gezeig, dass sich die Wasserstandsmuster im St.-Lawrence-Strom von Jahr zu Jahr wiederholen, aber auch Variabilität aufweisen. Die Methode war effektiv, um ähnliche Muster zu erkennen, selbst wenn sie in der Zeit verschoben oder in der Amplitude verändert waren. Für den Use Case der Hafen und Schiffahrtsplanung ist dies besonders relevant, da es ermöglicht, saisonale Anstiege im Wasserstand zu identifizieren und zu vergleichen. Eine Limitation würde dann entstehen, wenn extreme Wetterereignisse, langfristige Klimaveränderungen oder anpassungen am Flussverlauf wie zum Beipiel durch Dämme die Muster stark verändern. In solchen Fällen könnte die Cross-Korrelation welche hier verwendet wurde weniger zuverlässig sein und es müssten neue Messungen und Berechnungen durchgeführt werden.